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Abstract 

Many important properties of granular fluids can be represented by a system 
of hard spheres with inelastic collisions. Traditional methods of nonequilib- 
rium statistical mechanics are effective for analysis and description of the 
inelastic case as well. This is illustrated here for diffusion of an impurity 
particle in a fluid undergoing homogeneous cooling. An appropriate scaling 
of the Liouville equation is described such that the homogeneous cooling en- 
semble and associated time correlation functions map to those of a stationary 
state. In this form the familiar methods of linear response can be applied, 
leading to Green - Kubo and Einstein representations of diffusion in terms of 
the velocity and mean square displacement correlation functions. These cor- 
relation functions are evaluated approximately using a cumulant expansion 
and from kinetic theory, providing the diffusion coefficient as a function of 
the density and the restitution coefficients. Comparisons with results from 
molecular dynamics simulation are given in the following companion paper. 

PACS number(s):45.70.Mg, 45.70.-n, 05.70.Ln 
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I. INTRODUCTION 



It has long been recognized that rapid flow granular media have properties similar to those 
of ordinary fluids IHJ^. Attempts to sharpen this relationship have used idealized systems 
of hard spheres with inelastic collisions. Remarkably, the single feature of inelasticity allows 
reproduction of many qualitative features observed in real granular fluids. This observation 
is based on the derivation of hydrodynamic equations from kinetic theory 1^-0] , direct Monte 
Carlo simulation of kinetic equations 0,^, molecular dynamics simulation of dense fluids 
T0| , and controlled experiments on inelastic hard spheres [pA| -[T6|. Consequently, there 



is a growing opinion that the traditional methods of nonequilibrium statistical mechanics 
applied to such model systems provide the means to understand granular media at the most 
fundamental microscopic level. The objective here is to give further support for this opinion 
by a detailed application and test of linear response methods applied to diffusion of an 
impurity particle in a fluid of inelastic hard spheres. Some adaptation is required since the 
reference states are inherently nonequilibrium, but the central ideas of linear response for 
normal fluids are retained. A preliminary report of some of these results presented here has 
been given in . 



Diffusion is the prototype transport process and the associated diffusion equation is the 
prototype hydrodynamic description for macroscopic dynamics. For normal fluids, diffusion 
in a system of hard elastic spheres also has been the testing ground for many body methods 
in nonequilibrium statistical mechanics (density dependence, correlated many-body colli- 
sions, mode coupling, percolation, glass transition). The benchmarks have been set by 



accurate molecular dynamics simulation for normal fluids |18| and fluid-like lattice gas cel- 



lular automata |T^. In this way, the conditions for macroscopic diffusion and the accuracy 



of methods for predicting the diffusion coefficient are known over a wide range of densities. 
There is evidence based on kinetic theory and molecular dynamics simulation that similar 
studies are relevant for the system of inelastic hard spheres [^0|-p2| . On the other hand, there 
are significant differences to confront. An isolated fluid of inelastic hard spheres does not 
support a Gibbs state or any other stationary state, since the collisions lead to a continual 
loss of energy or "cooling" . Instead, the analogue of the equilibrium state is a "homogeneous 
cooling state" (HCS) whose time dependence is entirely given via the mean square kinetic 



energy of the particles |2^J2^J^,^. As shown below, the scaling property associated with 
this state allows a change of variables in terms of which a stationary, but non-Gibbs, state 
results. 

In the next section, the Liouville dynamics for inelastic hard spheres is reviewed 



T7| , |2^P^ . The corresponding nonequilibrium statistical mechanics is given in terms of 



the Liouville equation for the ensemble, and also in terms of the BBGKY hierarchy for 
the associated reduced distribution functions. Next, it is shown that the Liouville equation 
supports a scaling solution describing the HCS. A time dependent temperature is defined in 
terms of the mean square velocity, in the same way as for elastic collisions, hence the ter- 
minology "homogeneous cooling state" . A transformation to dimensionless variables allows 
a representation of expectation values and time correlation functions in terms of stationary 
state averages, just as for the Gibbs state. 

The probability density for the position of a tagged or impurity particle as a function 
of time is considered in Sec. |T|. Linear response methods (now inherently nonequilibrium 
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linear response) are applied to obtain a diffusion-like equation for long wavelengths, with 
a time dependent diffusion coefficient. The analysis parallels closely that for ffuids with 
elastic collisions, except that the dimensionless time (the collision number) required to ac- 
commodate the cooling is logarithmically stretched compared to real time. Consequently, 
the mean square displacement approaches linearity in the collision number rather than real 
time. The diffusion coefficient is expressed as a time integral of the velocity autocorrelation 
function which becomes a Green-Kubo expression for long times. This limit also establishes 
the time scale for the validity of the diffusion equation, and represents a clear example of 
macroscopic hydrodynamics for a granular system. 

An approximate evaluation of the velocity correlation function of the tagged particle is 
carried out in Sec. |V| by two different methods. The first is based on a leading cumulant 
expansion, while the second is an evaluation by means of kinetic theory. The two results are 
quite similar, and their relationship is clearly identified. Also, it is confirmed in Appendix 
that the evaluation of the Green-Kubo expression by kinetic theory agrees with that 
obtained by the Chapman- Enskog solution to the Enskog kinetic equation for the distribution 
function. Moreover, for elastic collisions these approximations are known to give an accurate 
description of the diffusion coefficient over a wide range of densities. This accuracy cannot be 
assumed a priori to apply for inelastic collisions, since the stationary state is nonequilibrium. 
The diffusion coefficient is given as a function of the density, restitution coefficient, and ratio 
of temperatures for the ffuid and impurity particles. This latter parameter is a peculiarity 
of granular ffuids for mixtures, where the lack of detailed balance leads to a HCS with all 
species having the same cooling rate but different temperatures [^. The consequences 
for impurity diffusion and mobility have been discussed elsewhere recently [^,^. The last 
section contains a summary and discussion of the main results. The theoretical developments 
presented here are tested and extended in a companion paper by comparison with molecular 
dynamics simulation results over a wide range of densities and inelasticies for the particular 
case of self diffusion 



II. STATISTICAL MECHANICS AND THE HOMOGENEOUS COOLING STATE 

The system considered is composed of a fluid of N identical hard disks or spheres (mass m, 
diameter a, and ffuid- ffuid coefficient of normal restitution a), and an additional impurity 
particle (mass mo, diameter ctq, and ffuid-impurity coefficient of normal restitution ao). 
The position and velocity coordinates of the fluid particles will be denoted by {qj,Vj;2 = 
1,...A^}, while those for the impurity particle are qo,vo. The dynamics consists of free 
streaming until a given pair of ffuid particles, i,j, or a fluid and impurity pair, 0,i, is in 
contact. At that point, the relative velocity of the pair changes instantaneously according 
to the inelastic collision rules 

gij = gij - (1 + «) (o- ■ gy)5-, goi = goi - (1 + ao) (^ ■ goi) ^- (1) 

Here gij = Vj — Vj and goj = vq — Vj are the relative velocities, and is a unit vector directed 
from the center of particle j to the center of particle i {i to 0, respectively) through the point 
of contact. The coefficients of restitution have values in the range < a,aQ < 1, measuring 
the degree of inelasticity. The special case of elastic particles is given by a = ao = 1. The 
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center of mass velocity is unchanged so that the total mass and momentum of the pairs are 
conserved in such collisions. However, there is an energy loss for each fluid-fluid particle 
collision 



TTl 

E,,-E,, = --il-a')ia-g,,f (2) 



and for each fluid-impurity particle collision 



Eo^ - Eo^ = -^{1 - al){a ■ goi)\ (3) 

where /x = mQ-m/ (mo + m). The state of the system at time t is completely characterized 
by the positions and velocities of all particles at that time and is represented by a point 
Tt = {qo(i^), • • • , (lN(t), vo(t), . . . , VAr(t)} in the associated 2d {N + 1) dimensional phase 
space, where d = 2 for hard disks and d = 3 for hard spheres. The sequence of free 
streaming and binary collisions determines uniquely the positions and velocities of the hard 
particles at time t for given initial conditions at t' < t. A more complete notation expressing 
this dependence on initial conditions is Tt(Tt'). Thus, just as in the case of elastic collisions, 
the microdynamics for this system corresponds to a deterministic trajectory in phase space. 

Observables of interest are represented by the same phase functions as for elastic colli- 
sions, A{rt), and their average for given statistical initial data at t = is defined by 

(A(t);0)= JdTp{T)A[T,{T)], (4) 

where p(r) is the probability density or ensemble for the initial state, normalized to unity. 
An equivalent representation of this average is obtained by changing variables to integrate 
over Tt rather than over F. This change of variables is possible since trajectories in phase 
space do not cross, and T can be expressed in terms of Tt denoted by T^^ (Tt). This allows 
the time dependence in Eq. (||) to be expressed in terms of the probability density 

{A{t);0) = JdTp{T,t)AiT)^{A;t), (5) 

with the probability density at time t given by 

p(r,t) = j[r-i(r),r]p[r-i(r)], (6) 

J (r, Tt) being the Jacobian of the transformation. 



A. Liouville Dynamics 



For practical purposes it is useful to identify the generators L and L for the two above 
representations, defined by 



{A{t);0)= I dT p{T)e'^A{T) = I dT e'^Xr) A{T) 



(7) 



The last equality is consistent with the adjoint relationship implied by Eqs. and (^. 
These are not the usual generators of Hamilton's equations for continuous forces, but are 
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somewhat more complex due to the singular nature of hard particles. Such generators have 
been discussed in detail for the case of elastic collisions and the analysis extends quite 
naturally to the inelastic case as well [0,^,^,^, with the results 



N 



N 



L = + vo ■ Vo + E^(0' 0, L = Lf + vo-Vo-Y. ^{0, i). 



i=l 



1=1 



(8) 



Here Lf and Lf are the generators for the fluid particles alone 



N 



1 



N N 



N 



N N 



1=1 i=l jj^i 1=1 i=l j^i 



(9) 



The terms involving spatial gradients generate free streaming while the others describe 
velocity changes. The binary collision operators T{i,j) and T{i,j) for particles i and j are 



cr 



dQQ{ 



S-)\gij ■ o-\6{qij - (T){hij - I] 



(10) 



T(z, j) = a^-i j dn e(g,, ■ a-)|g,, • a\[a-^5iq.,, - a)b-/ - 5(q,, + a)], (11) 

where dQ denotes the solid angle element for the unit vector a, cr = aa, is the relative 
position vector of the two particles, 6 is the Heaviside step function, and is a substitution 
operator, bijF{gij) = F(gjj), which changes the relative velocity gij into its scattered value 
gij, given by Eq. (|l|). On the other hand, it does not change the velocity of the center 
of mass of the two particles. The operator is the inverse of bij and characterizes the 
"restituting" collision. The binary operators for collisions between fluid particles and the 
impurity are similar to those for collisions among fluid particles, 

T(0, t) = a"-' J dfi e(-go. ■ a-) |go. • ^\6{qo^ - ^) {bo^ - 1) , (12) 



T(0, z) = a^-' I dneiSio ■ 5-)|go. ■ ^\[ao'5{qo^ - - 5(qo^ + ^)], (13) 

with = (cr + o"o) /2 and W = aa. The explicit forms for and 

1 + a 



= Sij ■ gij) (14) 

a 



^o/goi = got - (5- ■ goi) a. (15) 

The dynamics for the phase functions and the equivalent Liouville equation distribution 
function, p{T,t), follow directly from Eq. (|^): 

{dt~L)A{T,t) = (}, (16) 
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{dt + L) p(r,t) = o. 



(17) 



The BBGKY hierarchy of equations for the reduced distribution functions is obtained 
by partial integration of the Liouville equation over the position and velocities of N — I fluid 
particles, 



dt + L{xo, 



f 



i=0 



'^'^'\xo,...,xi+,-t), (18) 

with the reduced distribution functions defined by 

f+^\xo,...,xut)=N'j dxi+i...dxNp{{xi};t), (19) 

where Xi = {q^, Vj} denotes the position and velocity for particle i and L{xo, . . . ,xi) is the 
Liouville operator for a system of / fluid particles and the impurity particle. Moreover, the 
limit of large N has been considered. The above results in this section provide the basic 



tools and definitions of nonequilibrium statistical mechanics for granular media . 



B. Homogeneous Cooling State (HCS) 

Stationary solutions to the Liouville equation ( p!7D are expected when suitable external 
forces or boundary conditions are imposed. However, there is no stationary solution for an 
isolated system, corresponding to the spatially homogeneous Gibbs state, due to the inherent 
time dependence following from loss of energy in collisions. This can be seen by calculating 
the rate of change of the mean square velocity of a fluid particle in an isolated state. For 
purposes below, the latter is used to define a kinetic temperature according to 

T(t) = l(mv2;t) = imi;2(t). (20) 

In addition to T{t), Eq. (|^) defines the associated thermal velocity v{t) (a factor of Boltz- 
mann's constant usual in elastic systems has been deleted since there is no zeroth law of 
thermodynamics for granular media; alternatively, to incorportate the elastic limit a = 1 the 
temperature should understood as defined in units such that ks = 1). The time dependence 
of these quantities can be calculated using the Liouville dynamics to get 

dtT{t) = -mm, (21) 

where ({t) is the "cooling" rate due to inelastic collisions, 

C (t) = (1 - «')^^^^^ / / d^i J J dn e(gi2-a-)(gi2-a-)V^'^(qi, vi, qi + a, V2; t), 

(22) 

/*^^)(qi, vi, qi + cr, V2; t) being the reduced two fluid particle distribution function at contact. 
The latter is in general defined by 
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/(2)(a;i,X2;t) = I dxof^^\xo,xi,X2;t). (23) 

For an homogeneous system, its spatial dependence occurs only through qi2. Upon deriving 
Eq. (p^ we have taken into account that the contribution from the impurity particle is 
negligible in the limit of large N. 

In place of the Gibbs distribution, it is assumed that there is a homogeneous scaling 
solution phcs(X,t) to the Liouville equation, for which all time dependence occurs through a 
scaling of the velocity ("cooling") with the thermal velocity v{t) p3|J 



PHcsiT, t) = [ev{tr'^''+'^ ({%/^, ^Mt)}) . (24) 

The dimensionless distribution function p^^^ is invariant under space translations, with the 
coordinates scaled (arbitrarily) relative to i = {na'^~^)~^, which is proportional to the mean 
free path (n is the number density of particles). Therefore, phcs(X,t) represents a spatially 
homogeneous fluid. Substitution of Eq. (|24D into the Liouville equation gives 



1 ^ d 



The self-consistent solution to the coupled set of equations ( plf ) and (^) determines the 
homogeneous cooling state (HCS). It is the analogue of the Gibbs state for elastic collisions 
and reduces to it for a = 1. For a < 1, the exact solution is not known (it is not simply a 
Gaussian in the velocities as for the Gibbs state) but its existence is supported by results 
from Monte Carlo ^ and molecular dynamics simulations . 



The notation in Eq. ( plf ) does not make explicit the dependence of phcs on parameters 
specific to the impurity particle. Since all velocities have been scaled relative to the thermal 
velocity determined by the fluid, there is a clear dependence on the mass ratio m/mo, 
expected from equipartition of energy. However, equipartition is strictly a property of the 
equilibrium Gibbs state and does not apply for the HCS Consequently, in addition to 



a dependence on the mass ratio, there is also a dependence on the temperature ratio T/Tq, 
where Tq is the temperature parameter for the impurity. The detailed form of the relationship 
between To(t) and T{t) will be discussed later on, but it follows from the condition that the 
cooling rates of the fluid particles and the impurity particles are the same in the HCS, 
implying that T /Tq is time independent. Through this section, the dependence of phcs 
on time independent parameters of the impurity particle will continue to be suppressed, 
although it will become important in the subsequent discussion of impurity diffusion. 

Some interesting consequences follow from the velocity scaling of the distribution function 
associated to the HCS. The reduced distribution functions also have this property, so it is 
easily verified from Eq. (^ that C {t) oc T^/^{t). Then, Eq. (|^) can be integrated for the 
explicit time dependence of T(t), 



T{t) = T{t') i+\am-t') 



(26) 



The temperature is seen to have an algebraic decay in real time (Haff's law ||T|). For the 
analysis of the HCS, it is more convenient to use the dimensionless time scale 
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s{t,t') = / dTu.ir), 
Jt' 

where i^dt) is an average collision frequency given by 

u,{t)=v{t) /i. 



(27) 



(28) 



Thus s{t,t') is a measure of the average number of collisions per particle in the interval 
{t',t). The integral in Eq. (|27|) can be performed using Eq. (^Bj) with the result: 



In 



1 + it - f) 



where we have introduced the dimensionless cooling rate 

C = iC{t)/v{t). 



(29) 



(30) 



It follows from dimensional analysisis that (* is time independent. The cooling in terms of 
the dimensionless time s is exponential, 



T{t) = T{t')e 



-Cs{t,t') 



and, consequently. 



(31) 



(32) 



Since C is proportional to (1 — a^), there is a crossover from logarithmic to linear relationship 
between the two time scales for weak inelasticity. 

Knowledge of the time dependence of v{t) also implies that for many average properties. 
For example, the average value of A{T) in the HCS can be written 



{A-t) 



hcs 



drpU^,t)Air) = J drpUr)A{{eq:,v{t)^r*}). 



(33) 



Use has been made of the scaling form (p^ and T* = {q^, v*} = {qi/i,Vi/v{t)}. This last 
result suggests that the transformation to dimensionless form may admit a stationary state 
representation for the HCS. To see that this is the case, define for a general state p(T,t) 



p{T,t) = [iv{t)]-'^''^''^p* {T\s) 



(34) 



where we have introduced the same scaling for space, time, and velocity as above. Substi- 
tution of this into Eq. (|l^ gives the dimensionless Liouville equation 



with the definitions 



+ p*{r,s) = o, 



C* = lc [IC* + d{N + 1)] + L* 



(35) 



(36) 



' L 



{qi=q*,Vi=v*} 



(37) 
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TV Q 

^* = Ev;-^. (38) 

This transformation of the Liouville equation exphcitly accounts for the colhsional coohng, 
and in this form stationary solutions are now possible. In fact (|25|) becomes 

^*pL = 0, (39) 

so the dimensionless HCS is a stationary solution to (^). 

The average value of A{T) for a general state p(r,t) becomes 

{A-t) = J dTp{T,t)A{T) = J dT*p*{T*,s)A{{iqlv{t)^r*}) 



dT* 



^ J drp*{T*, 0)U{s, 0)A{{icC, v{t)v*}), (40) 

where U{s, s') and U{s, s') obey the equations 

{ds + jC*)uis,s') = 0, {ds- C*)Uis,s')=0, (41) 

with the initial conditions U{s', s') = U{s', s') = 1. The new generator for the dynamics of 
the phase functions is 

^* = ^C*/C* + L*, ^* = [i^]{,..,*,.....|. (42) 
For the special case p*(r*,0) = pL.(r*), Eq. (H) reduces to (H). 



The stationary representation is the most natural one for both theoretical developments 
and for computer simulation, as is discussed in the following companion paper Similarly, 
the physically relevant time scales are those expressed in terms of the average collision 
number s rather than the real time t. It is appropriate at this point to note that although 
p^^g is a stationary solution to (§^), there is convincing evidence from both theory and 
simulation that it is unstable to long wavelength spatial perturbations and spontaneous 
fluctuations P^ , P3[] . In the following sections, time correlation functions are considered for 
the HCS and use is made of stationarity and spatial homogeneity. The results must be 
understood as applying to system sizes for which the instability does not occur, or on time 
scales that are short compared to those required for growth of spatial structures. 



C. HCS Averages and Correlation Functions 

The HCS time correlation function for two phase functions A{r) and -B(r) is defined as 

CABit,t') ^ {Ait)Bit'y,0) - (A(t);0)(fi(O;0), (43) 

with t > t' > 0. Here and below the brackets (; t) denote an average over the HCS at time 
t. For a system with elastic collisions in equilibrium, the above autocorrelation function 
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can be reduced to a single time correlation function, using time translational invariance and 
stationarity of the Gibbs state. In the case of inelastic particles, the HCS is not stationary, 
but the scaling property p^ ) can be used to transform the correlation function to an effective 
time stationary average. First, use time translational invariance to write 



{A{t)B{t'y,0) = {Ait-t')BiOy,t'). 
Next, transform to dimensionless variables to get 



(44) 



{A{t)B{t'y,0) = J dTphU^,t')A{t-t')B{0) 



^(Uq:,^(Ov*}) 5({^q-,^(0v:}), (45) 



where (|33| ) has been used. Next note the identity 



e 2^ 



The correlation function now can be written 

{A{t)B{t'yo) = J rfrvL.(r*)[f/(t,t')^({H*,^Wva)]i?(H*,^(Ov*}), 

where 

U{t,t') = e'^<=(*')(i-i')i*e|C*s(</)/c*^ 



(46) 



(47) 



(48) 



This time evolution operator can be identified by differentiating with respect to s{t,t'), 
taking into account that in the HCS (* is time independent. 



dU{t,t') 
ds 



Uit,t') 
U{t,t') 



ds ^ ' 2^ 



Consequently, 



U{t,t')C 



(49) 



(50) 



Note that this is propagator is the same as that in Eq. ^ specialized to the HCS. In this 
case C* becomes time independent, allowing the simple exponential representation. The 
correlation function now can be written in the final form 



(A(t)i?(t');0)= /rfrpL.(r 



(51) 



This is a primary result of this section. The time correlation functions depend on the 
dynamics through the collision number s{t,t'). All additional time dependence occurs triv- 
ially through the thermal velocity. This is most evident when A and B are homogeneous 
functions of the velocity. 



10 



^({£q*,t;(t)va) = v^mHi^l^^:}), S(£q*,t;(Ova) = v\t')B{{iql^r:}). (52) 
Then the correlation function becomes 

{A{t)B{ty,0) = v^{t)v\t'){A{s)Br, (53) 

{A{s)Br = J rfrXc.r)A({£q*,va,3)i?(w,va), (54) 

and the phase function y4({£q*, v*}, s) is 

v*}, s) = e^(*-*')^* A({^q*, V*}) (55) 

This stationary state representation for the time correlation functions simplifies considerably 
the analysis of response functions in the next section. 



III. IMPURITY PARTICLE DIFFUSION 

In this section, the diffusion equation and associated expressions for the diffusion coeffi- 
cient are derived for a granular system in the HCS. The probability density P(r,t) to find 
the tagged or impurity particle at point r at time t, given it was at the origin at t = 0, is 
defined by 

P(r, t) = V{6[qo{t) - r]5(qo); 0) = (5[qo(t) - qo - r]; 0), (56) 

where the angular brackets indicate as above an average over an initial HCS and V is the 
volume (for d = 3) or surface (for d = 2) of the system. The second equality is a consequence 
of the translational invariance of the HCS. The conservation law for probability follows by 
differentiation of Eq. (^) with respect to t, 

atP(r,t) + V- J(r,t) =0, (57) 

with the probability flux J(r, t) identified as 

J(r, t) = (vo(t)5[qo(t) - qo - r]; 0). (58) 

The interest here is in the limiting behavior of Eq. (0) in the hydrodynamic regime, 
which corresponds to the long wavelength region. The long wavelength spatial dependence 
of J(r, t) can be obtained from a Fourier representation of P{r,t), 

P{k,t) = J dre'^''P{r,t). (59) 

To get a formal equation for P(k, t), it is useful to introduce the index of the distribution 
Cik,t) by 

P(k,t) =6^^^'*), (60) 

i.e.. 
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C(k, t) = ln(e^''-['i°(*)-'i°(o)l ; 0) . (61) 
Differentiation with respect to time of Eq. ( |60|) yields 



dt-C{k,t) 



P(k,t) = 0, (62) 



where the dot over C denotes the derivative with respect to time. For long wavelengths 
{ki <^ 1), C{k,t) can be expanded to order k"^, 

~ _ (ik- vo(t)e^'^-tqow-q»w];0) 

^-^(vo(t)-[qo(t)-qo(0)];0). (63) 

Substitution of this expression into Eq (^) and inverting the transform gives Eq. (0) with 
the identification of the probability flux as 

J(r,t) = -D(t)VP(r,t). (64) 

The flnite time diffusion coefficient is 

^W = ^|(|qo(t)-qo(o)r;0). (65) 

This will be referred to as the Einstein form, relating the diffusion coefficient to the mean 
square displacement. The equivalent Green-Kubo form, in terms of the velocity autocorre- 
lation function (VACF), is derived by using the relationship 

qo(t)-qo(0)= tdt'Mt'), (66) 
Jo 

with the result 

D{t) = ^ fdt'{Mt)-^oit');0). (67) 
a Jo 

For normal fluids with elastic collisions, the diffusion constant follows from the long time 
limit D = limt^oo D{t), or equivalently from the coefficient of the mean square displacement 
when it becomes linear in t. This limit occurs for times large compared to the mean free 
time. Since the latter is time dependent in the HCS, the usual conditions to establish a 
diffusion constant, and consequently the diffusion equation, must be modified for granular 
media. This is done by introducing the dimensionless diffusion coefficient 



Using the representation ( |55| ) for the correlation function in ( p7D gives the Green-Kubo form 
ClS db stationary average 
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D*{s) = t dt'v{t)v{t'){^r*{s) ■ v*y 



d&u, (t) Jo 



d Jo 



ds'i^oi^') ■ V*)*. (69) 



d 

1 rs 



d 

In going from the first line to the second hne we have written s{t, t') = s{t, 0) — s{t', 0) = s—s'. 
Similarly, from Eqs. (|65D the corresponding Einstein form is 

These are the stationary average representations for the diffusion coefficient, and are the 
primary results of this section. In terms of the time scale s, the mean square displacement is 
expected to become hnear for s 1, and the velocity autocorrelation function is expected to 
decay to zero also on this time scale. The physical interpretation of this limit is the same as 
for elastic collisions, since s is essentially the number of collisions per particle. However, due 
to the time dependence of the collision frequency, the correlation functions are expected to 
have the proper behavior with respect to s rather than t. This will be shown more explicitly 
in the next section. The dimensionless form of (|57|) and ( |64D reads 

dsP*{r*, s) - D*{s)V*^P*{r*, s) = 0, (71) 

where P*{r*, s) = i'^P{r, t). Clearly, this becomes the usual diffusion equation for sufficiently 
large s, if D*{s) tends to a constant. 

It is useful to introduce a dimensionless VACF for the impurity particle that is normalized 
to unity at s = 0. This requires calculation of (wq^)*. It might appear from Eq. (^Of ) that 
this can be obtained simply in terms of the mass of the impurity and the temperature of 
the fluid. However, it has been shown elsewhere that mechanically different particles in a 



common HCS do not have the same temperatures |^6[, as already mentioned below Eq. 
Thus < fg^ >* is given by Eq. (^) with both m and T{t) replaced by mo and To(t), 

^9 ^ dTnit) d , 

where (phcs is the ratio of the square of the thermal velocity for the impurity particle relative 
to that for the fluid particles, 

'^^^ = ^^^^ 

Since the cooling rates of the fluid, Cif), and the impurity particle, Co(^)5 are the same 
and they are proportional to the square root of the temperature, the above ratio is time 
independent. The condition for equal cooling rates also determines (t)hcs- These rates are 
calculated to good approximation using a local equilibrium ensemble in Appendix 0, and 
are given by 
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c-Wm^^'-""'' '''' 

Q = ^- (l - hl±^] (1 + (75) 

\ <Phcs ) 

where /i = (1 + a^) m/2 (m + mo), and v* is a dimensionless impurity particle collision rate, 

The factors x ^'^^ Xo ^"^^ the pair correlation function for the fluid-fluid and the fluid- 
impurity particles at contact, respectively. A more accurate calculation for the case of hard 
spheres {d = 3) is carried out in Ref . . Equating Eqs. (|7^ and (|75D provides the equation 

for (f)hcs 

{l + <f>.csf'(l-h^^] = ^-^. (77) 

V <Phcs J 

This gives a cubic equation which has a unique real, positive solution for all allowed values of 
h and C* 1^*- Foi' elastic collisions, (phcs rn/rriQ as required by the equipartition theorem. 
Qualitative changes in this solution, similar to a phase transition occur in the limit /z — ^ 
p8| , but will not be discussed here. The normalized VACF is now given by 



CUs) - ^^l^ = is) ■ <Y (78) 



and the diffusion coefficient in Eq. (p9| ) becomes 



s 



D*is) = ^ I ds'Clis'). (79) 







IV. APPROXIMATIONS 



In the following, two approximations, originally developed for fluids with elastic collisions 
31], are applied to calculate the VACF in the HCS. The first method uses a leading order 



truncation of a cumulant expansion, while the second uses an approximate kinetic equation. 
The results confirm the expected time scale for transition to hydrodynamics and provide 
the detailed dependence on density and restitution coefficients. 



A. Cumulant expansion 

The cumulant expansion of the VACF is 

OO 



Cl^{s) = exp 



p=i p- 



(80) 
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The coefficients 00* are determined from the initial time derivatives of the correlation func- 
tion. Clearly, truncation of this expansion at any order is asymptotically exact at short 
times and also for small h (heavy impurity) since each time derivative contributes a factor 
of h. The simplest such approximation retains only the leading term 



c: 



(81) 



with 



d_ 

ds 



InC^^(s) 



< (^*v* 



s=0 



^hcs 



(82) 



Use has been made of the definition of the s dependence in (^8]). The corresponding ap 



proximation for the mean square displacement is obtained by integrating Eq. (81) twice 



qo 



■'hcs 



1 



(83) 



For elastic collisions, this approximation is known to be accurate for short as well as long 
times, and for a wide range of densities and mass ratios. The resulting diffusion coefficient 
and diffusion constant are then found by substituting Eq. (|81|) into Eq. (|7£ 



D*{s) = 1 - e 



(84) 



D*=\imD*(s) = ^. (85) 

Clearly cj* is a characteristic dimensionless collision frequency for the impurity particle. The 
ffist cumulant approximation conffims the expectation that the mean square displacement 
becomes linear in s and the velocity autocorrelation function decays for s >> ul- Conse- 
quently, the macroscopic diffusion equation applies on this time scale as well. The collision 
frequency ujI is evaluated in Appendix ^ by using a local equilibrium ensemble, with the 
result: 

ujl = ~C + i^*{l + (^hcsf\ (86) 

where C*and z/* are defined in ( ffil) and (|76D , respectively. Substitution of this into Eq. (|85D 
gives 

D* = . (87) 

It is possible to show that D* is positive and finite for all values of the density and restitution 
coefficients. 
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B. Kinetic Theory 



Perhaps the most accurate and detailed evaluation of time correlation functions is via 
kinetic theory methods. These can be applied as well to the case of inelastic collisions 
Il30| , |35|jr?| . To show this, first use the adjoint property of the Liouville operators to write 
the velocity autocorrelation function in the form 



CL(s) 



^hcs 



^hcs 



Next, define the dimensionless reduced correlation functions in a manner similar to (|T 



X, 



0' 



so that the VACF can be written as 



■'hcs 



J dx*Q V* • V 



Xq, s). 



(88) 



(89) 



(90) 



It is easily verified by direct differentiation that the -j/''-'^^'' functions obey a hierarchy of 
equations similar to (IH), the first of which is 



{^, + ^r*.Vl)^^^'^'\x;,s) + ^C 



d 



2^ dv* 



dxlT''{0,lW^'^\xl,xl,s). (91) 



A kinetic equation results from a closure approximation in the above equation that 
expresses i/?**-^'' as a functional of Formally this is possible since both iI}*^^\xq, s) and 

if)*^'^\xQ, xl, s) are linear functionals of xI)*^^\xq, 0). In principle, this functional relationship 
can be inverted to give ip*^^\x*f^,xl, s) = **^^^ \xl,xl, s\il)*^^^]. Use of this in Eq. (|T]) 

provides the closed kinetic equation for In practice it is a difficult many body problem 

to discover this functional. However, its form is easily calculated at s = 0. Equation (|89D 
gives 

[x*,xi,s = o\r^'\o)] = r(2)(x*,xi)v* = r^'\x*,xi) [f<'\x*)Y'r^'\xio), (92) 

where f*^^\xQ) and f*^'^\xQ,xl) and are the reduced distribution functions associated with 
Phcsi^*)- -^"^^ fluids with elastic collisions, the approximation 



(93) 



is accurate over a wide range of low to moderate densities. The same approximation in the 
current case leads to the kinetic equation 



(a. + v*-vS)V*^^)(^S,^) + 



C* d 



dxlT*iO,l)f<'\x*,xl) \r^'\xl)r r^'\xl,s). (94) 
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This equation is exact at short times by construction and its use for longer times can be 



interpreted as a Markovian approximation [^. For the VACF only the spatial integral of 
iI:*^^\xq, s) is required, so the final representation of Eq. (pOD becomes 

C^:.(^) = ^/rfvSvS>*K,^), (95) 

where 

n<,s) = J dciir^'\xl,s) (96) 

obeys the kinetic equation 

ds^K, s) + ■ [v*V *(vS, s)] = irK, s). (97) 

The collision operator I is defined as 

Jr«s) ^ J dxlT\QA)f<'\xlxl) [r(i)(a;S)]"V(vS,3). (98) 
The kinetic equation ( P7| ) has to be solved with the initial condition 

V'*(vo,0) = rW(vo)vS, (99) 



as follows directly from Eq. (^). To obtain Eq. ( pTf ) from Eq. (0), use has been made of 
the fact that the HCS distributions are spatially homogeneous. For a = 1 the linear collision 
operator I is non-negative, and the correlation function is found to decay on a time scale 
of the order of s ~ the initial rate of decay. Although not proven, it is reasonable to 

assume that the spectrum of X is qualitatively similar for a <1. Since the kinetic equation 
is exact at s — 0, the leading term in a cumulant expansion of Eq. (|89D also is exact and 
agrees with Eq. (|HI|). More generally, this approximate kinetic equation gives contributions 
to all higher terms in the cumulant expansion. However, for a = 1 these corrections are of 
the order of a few percent except when the size or mass ratio of fluid and impurity particles 
differs greatly from one. 

The diffusion constant is obtained by integrating (|90|) to get 

^* = -^/rfv*v*-XK), (100) 

where X (vq) is the solution to the integral equation 

^X.. (vo) - ■ i^o^^ {<)] = r^'^ i<) <■ (101) 

The only (left) eigenfunction with vanishing eigenvalue of the operator defined by the left 
hand side of the above equation is 1. The right side is orthogonal to 1 so the Fredholm 
alternative (solubility condition) is satisfied and a solution to this equation exists. It is 
shown in Appendix |^ that the diffusion coefficient given by Eq. ( |100| ) is the same as that 
obtained from the Chapman-Enskog solution to the Enskog kinetic equation, if velocity 
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correlations are neglected for f*^'^^ (xq, x*) in the definition of X, Eq. (pSl). Recent computer 
simulations confirm that such neglect is a good approximation . Finally, if Eq. ( |101D is 
solved as an expansion in Sonine polynomials (the usual method for elastic collisions), the 
first approximation yields again the leading cumulant approximation. Thus the cumulant 
approximation, linear kinetic theory, and Chapman-Enskog solution to the Enskog kinetic 
equation are all the same modulo small differences due to velocity correlations and higher 
order Sonine polynomials. 



V. DISCUSSION 

In this paper, it has been shown that standard linear response theory can be extended in 
a natural way to describe diffusion in a system of inelastic hard spheres in the homogeneous 
cooling state. The response functions of the system are given in terms of stationary state 
averages corresponding to a an effective dimensionless dynamics. In particular, the dimen- 
sionless time scale is related to the average number of collisions per particle taking place in 
the system, and it is the physically relevant one to analyze the aging to a hydrodynamic 
stage. Similar Einstein and Green-Kubo expressions have been obtained for the shear vis- 
cosity |]38|j39| , which is essentially a diffusion process like that considered here. Response to 



an external force has been studied to determine expressions for the mobility |2^. Transport 
coefficients associated with longitudinal hydrodynamic modes (e.g. thermal conductivity) 
pose special problems and will be discussed elsewhere. Beyond their formal interest and util- 
ity for approximate analysis, the results derived here enable a systematic nonperturbative 
study of transport processes in granular fluids by means of molecular dynamics simulation 
of the response functions, just as for fluids with elastic collisions. 

It is interesting to consider the particular case of self-diffusion, i.e. when the impurity 
particle is mechanically equivalent to the fluid particles, for which previous analysis have 
been carried out. The expression for the self-diffusion coefficient in the first cumulant ap- 
proximation is obtained by considering the limit (j)hcs = 1, /i = (1 + a)/4, and ao = a in Eq. 
(13). This yields 

dT(d/2) , , 

(1 + a)2x2i/27r— 



As expected, this result agrees with the expression derived in Ref. ||2^ from a Chapman- 



Enskog solution to the Enskog-Lorentz equation in a first order polynomial expansion. A 
previous approach to self-diffusion in the HCS, also in the context of linear response theory, 
has been developed by Billiantov and Poschel 0]. Their result differs from Eq. ( p.02|) in 
a factor of (1 + a)/2. Quite peculiarly, the same result had been derived before by Hsiau 
and Hunt and Savage and Dai independently, from approximate solutions of the 
Enskog-Lorentz equation. The discrepancy between the result in Ref. and the one 
reported here is not due to a different degree of approximation, but has a fundamental 
physical origin. Brillantov and Poschel assume that the velocity autocorrelation function of 
the tagged particle has the form 



\C,,{t,i!)\^p = v{tye (103) 
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with the relaxation time tb being inversely proportional to the initial slope, in the actual 
time scale t, of the VACF. This is to be contrasted with our analysis, based in the relevance 
of the time scale scale defined by the number of collisions. 



Cyv{t,t') = v{t)v{t') exp 



(104) 



with a;* determined from the initial slope of the VACF in the dimensionless scale defined 
by Eq. (P?]). Making clear the crucial role played by this latter time scale for the study of 
response functions is one of the main goals in this paper. 



The quality of the simple approximations given here is studied in the following paper 
by comparison with molecular dynamics simulations for both the Einstein and Green-Kubo 
forms. Only the case of self- diffusion is considered. It is found that the agreement is very 
good for low densities and all degrees of dissipation, but there are large discrepancies at 
high density and large dissipation. The possible reasons for the discrepancies are discussed 
there. 



VI. ACKNOWLEDGMENTS 

The research of JWD was supported by National Science Foundation grant PHY 9722133. 
JJB acknowledges partial support from the Direccion General de Investigacion Cienti'fica 
y Tecnica (Spain) through Grant No. PB98-1124. JFL acknowledges support from the 
Universite Libre de Bruxelles. 



APPENDIX A: COOLING RATES IN THE HCS 

In this Appendix, the cooling rates for the fluid and the impurity particle, as well as 
first cumulant for the VACF, are calculated from the Liouville dynamics of the system. For 
simplicity, a local equilibrium ensemble approximation, known to lead to accurate results 
1^, will be considered. The dimensionless cooling rates for the fluid and impurity particle 
are defined as 

^ v{t) dt r ^ ' ' ^ ^ 



Co 



dlnTn 



■(L*v*2)* 



(A2) 



v{t) dt diphcs 

respectively. For the derivation of the last equalities, use has been made of the property 



that in the HCS, for any dynamical variable AiV) having the scaling property (52) it is 

(A3) 



{LA-t) = ^^{L*A{{lcC,,^:W. 



The first cumulant for the VACF is given by Eq. (p2|) or, equivalently, 
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C 2 

In terms of the dimensionless distribution functions 

f<^){xixi) = N j dx;dx;...dx%pi,,ir), 

and 

X2) = A^^ J dx^dx^ . . . dx*]ypX^g{T*), 

the above quantities can be written as 

2 



Co 



dN 



■'hcs 



J dx\dxif<'^ (xt,x;)r*(i,2)vf, 



dxldxlf<'^ {xlx\)T\Q,l)wf 



(A4) 
(A5) 
(A6) 

(A7) 
(A8) 



^hcs 



dxidxif<^^ (x*,xDv*-[r*(o,iK]. 



(A9) 



The action of the binary coUision operators can be performed for further simpUfication, 
leading to 



1 



2dn* 



1 „*d~l 
I — 01 cr 



j dvt j dwl I rffi/*(2)(v*,V2*,q*2 = o-*)e(-f 



;i2 ■ '^j|gi2 • 



(AlO) 



4h 



■'hcs 



X 



2G*i ■ a + 



m + mo 



1 - ao) ■ 



(All) 



wi = — - 



4/i 



^hcs 



X ( G*i ■ a- + 



mo + m' 



(A12) 



where h is defined below Eq. (175|), is the reduced volume or surface of the system, and 
Ggi is the reduced velocity of the center of mass for particles and 1, i.e. 



01 



moVQ + m\l 

mo + m 



(A13) 



The above results are still exact. To compute the integrals two approximations are intro- 
duced. First, the velocity correlations in /^^^*(v^, Vg, = c*) and f*^'^^ (vg, v]',qoi = W*) 
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are neglected. Note that this approximation is required only for particles at contact and 
for the precoUision hemisphere. Significant velocity correlations exist on the opposite hemi- 
sphere and at larger separation of the particles, and no restriction is placed on these con- 
figurations. The second approximation is to represent the one particle distributions by 
Maxwellians. This is known to be a good approximation except for conditions of large me- 
chanical differences between the fluid and impurity particles. These two approximations are 
equivalent to write 

r(2)K,vt,q*, = 0^) =n*Xo(^)v^SK)^*K), (A14) 



r(2) (vt, v;,qt2 = ct)= n*\{a)v*{^rl)^*{^r;), (A15) 

where Xo(^) ^^nd xi'^) ^''^^ fhe fluid- impurity and and fiuid-fiuid pair correlation functions 
for particles at contact, and 

VoiO = ^e"^, (A16) 



IT 



d/2 



(A17) 



Now the integrations in Eqs. ( A10 )-( |Al2| ) can be performed. The calculations are straight- 
forward but lengthy, and will be not reproduce here. The results are given by Eqs. (|7^ -(|76D, 
and dH). 



APPENDIX B: ENSKOG KINETIC THEORY 



In the main text sections, linear response has been applied to describe diffusion directly 
from the Liouville dynamics. It is useful to see that equivalent results, in appropriate 
approximations, follow from kinetic theory as well. This is illustrated in this Appendix using 
the Enskog kinetic theory, expected to be valid over a wide range of densities and restitution 
coefficients. To review its origin and applicability to granular fluids, the first hierarchy 
equations of ([18|) for the impurity particle reduced distribution function /o^''(qo, vq, t) and 
one particle fluid reduced distribution function /'^^^(qi, vi, t), are written exphcitly, 

{^t + ^ro■Vo)Jf\qo,^r,,t)=a^-' Jdq, j dv^ j rffi e(-goi • ^) |gor^| 

X ao^5(qoi -^)&oi^ - 5(qoi + o^)] /^^^qo, vo,qi, vi,t), (Bl) 



;i2 ■ crj|gi2 • cr| 



a '(5(qi2-o-)6i2'-5(qi2 + <T)]/(2)(qi,vi,q2,V2,t), (B2) 

were /'^^^(qo, Vq, qi, Vi, t) and /'-^•'(qi, Vi, q2, V2, t) are the two particle reduced distribution 
functions for the impurity and one fluid particle, and for two fluid particles, respectively. 
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Their definitions are given by Eqs. ([T9|) and (0). A closure of tlie liierarcfiy is obtained by 
replacing in the above 

/^^^qo, vo,qi, vi,t) xo(qo,qi;i)/o^^(qo, vo,t)/^^nqi,vi,it), (B3) 

(qi, vi, q2, V2, t) -> x(qi, q2; t)/^'nqi, Vi, t)/^^) (qs, v^, t). (B4) 



The approximation is a generalization of that in and also of that in Eqs. ( [A14| ) and 



( |A15|) . It can be understood in two different ways. It is the exact Markovian limit if the 



initial distribution functions have the forms (p3|) and (B4) on the precoUisional hemisphere 
(exact at short times). It also follows if velocity correlations on this hemisphere are destroyed 
between collisions (Boltzmann's argument). In any case, it is known to provide a good 
description of the hard sphere fluid over a wide range of densities and times for elastic 
collisions, and this validity appears to hold as well for inelastic collisions. As in the elastic 
collisions case, the functions Xolqc^qiS^) and x(qi!q2;^) are taken to be the equilibrium 
configurational pair correlation functions as functionals of the nonequilibrium density. 

The Enskog approximation converts Eqs. ( |B^ ) and ( |B4D into a pair of kinetic equations. 
The kinetic equation for the fluid is autonomous while that for the impurity distribution is 
a functional of /'^^^(qi, vi, t), 

{dt + vi ■ Vi) /(qi, vi, t) = Je [qi, Vi|/(t)] , (B5) 

(9i + vo ■ Vo) /o(qo, vo,t) = Jel [qo, vo|/o(t), /(t)] . (B6) 

Here and below the superscript on the single particle distributions will be suppressed for 
simplicity of notation. The Enskog and Enskog-Lorentz collision operators are given by 

^i?[qi,vi|/(t)] ^a'^-^l J dQQ{-gu-^)\gi2-^\ 

X {a^^X [qi, qi - (^\n{t)] &r2V(qi, vi, t)/(qi - cr, V2, t) 

-X [qi, qi + cr|^W] /(qi, vi, t)/(qi + a, V2, t)} , (B7) 

JEL[qo, vo|/o(t), fit)] = cf-' J rfvi J dQ e(-goi ■ ^)|goi ■ 5-| 

X {"0^X0 [qo, qo - '^\nit)] 6oiVo(qo, Vo, t)/(qo - o", Vi, t) 
-Xo [qo, qo + '^\nit)] /o(qo, Vq, t)/(qo + ct, Vi, t)} , (B8) 

respectively. The HCS is a special case for which Eqs. ( |B5p and (|B^ ) reduce to 

dtfhcs = JE[Vl\fhcs{t)], (B9) 

dtfo,hcs = Jel [vo|/o,hcs(t), fhcsit)] , (BIO) 

with the collision terms 
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M^ilfhcsit)] = [ dv2 f dne{-gi2 ■ 5-)|gi2 ■ o- [(a - l) /hc.(vi,t,)Ac.(v2,t 



(Bll) 



JEL[vo\fo,hcsit), fhcsit)] =0^-^X0 J dVi J rff]e(-goi ■ CT)|g01 " ct| 

X («o '^w - 1) /o,hcs(vo, t)/hc.(vi, t). (B12) 

Use has been made of the fact that x [qi, qi — cr|n] = x and Xo [Qo? Qo ~ "'l'^] = Xo are 
independent of space coordinates due to translational invariance. Furthermore, the reduced 
distribution functions fhcs and foMcs have the scahng forms, inherited from (^), 

/o,hcs(vo,t) = [^t^(t)]-'/o%cs [vo/t^(t)] . (B13) 
Substitution of these expressions into the kinetic equations leads to 



d 

1 d 

2^°9v5 ' ^ vol/o*/,^^,/^^, . (B14) 

where the scaled position and velocity variables defined below Eq. (|33|) have been again 
introduced. The solution to these equations has been discussed elsewhere and will not 
be considered further here. 

Now consider more general spatially inhomogeneous states for the impurity particle, with 
the fluid still in the HCS. The probability density to find the impurity at a point qo in terms 
of the reduced distribution function is given by 

P(qo,t) = J rfvo/o(qo,Vo,t). (B15) 

Integration of Eq. ( P6|) over the velocity vq gives the conservation law for probability 

aiP(qo,t)+Vo- J(qo,t) =0, (B16) 



J(qo,i) = j c?VoVo/o(qo, Vo,t) (B17) 

which, of course, are the same as Eqs. (|57|) and (|58|) - The diffusion equation is obtained 
from a "normal" solution to the Boltzmann-Lorentz equation, where all space and time 
dependence occurs through P(qo,t) and T{t). The linear relationship ( |B15|) implies that 
such a solution has the form 

/o(qo,Vo,t) = P(qo,t)/i[vo/t;(t)]. (B18) 

The Chapman-Enskog method represents a normal solution to the kinetic equation as an 
expansion in the gradients. 
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/o(qo, vo, t) = fi"^ (qo, Vq, t) + e/^^^ (qo, Vq, t) + . . . , (B19) 

where e is a formal parameter representing the order of the spatial gradient. Similarly, the 
time derivative is obtained as an expansion in the gradients via the conservation equation, 

dt = 9f ^ + eSf ^ + . . . (B20) 

and the kinetic equation is written 

{dt + evo ■ Vo) /o(qo, vq, t) = Jel [qo, vo|/o(t), fhcsit)] . (B21) 

Substitution of Eqs. ( |B19| ) and ( [B20| ) into Eq. ( |B21| ) gives to zeroth order in the spatial 
gradients 



(0) .(0) 







J. 



EL 



^.0,^Q\fo\t)Jhcs{t) 



which has the solution 



/o^°^(qo,vo,t) = KP(qo,t)Ac.(vo,t). 



(B22) 



(B23) 



This gives J'-'^-' = and, consequently, dt^^ P = 0. The first order correction f^^^ is determined 



from 



5fVo^'^ + id't" + vo • Vo) = Jel qo, Vo|/o^^^(t), 



.(0) 



(B24) 



The contribution from df^^ vanishes, since it is proportional to d^^^ P{q^o,t) = — Vo • 
j(°''(qo,t) = 0. Also, the time derivative d^^^ can be expressed in terms of the cooling 
rate (. Using once again the dimensionless variables defined in the main text, Eq. ( P24|) 
becomes 

^* ^ - /oV.Vo ■ VSP*(qS,t), (B25) 



-^-ELJo - 



2 



vS/o 



where Tel is the linear Enskog-Lorentz collision operator, 

X^^/*« ^Jd^rlJ dQQ{-g;, ■ a)\g;, ■ a\ {a^%,' - l) rK)/o^'^(qS, vS,t), (B26) 



that is the same as the collision operator X in Eq. (pSf ) if velocity correlations, associated to 
the two particle distribution, are neglected in the latter. The solution to Eq. ( |B25| ) can be 
written in the form 

/o^'^(vo,q*,t) = XK)-V*P*(q*,.), 
where X(vq) is the solution to the integral equation 

C* d 



IelX, (v*) 



2 d^l 

Finally, the diffusion coefficient is identified from 

J^'Hqo,t) = / rfvoVo/r(qo, Vo,t) = -DVoP(qo,t) 



(B27) 
(B28) 



to get 



D* 



■'hcs 



(B29) 
(B30) 



These are the same results as those from Sec. Eqs. (|100|) and (|101|) . 
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